IFJPAN-IV-2006-7 

Propagation of Uncertainty in a Parton Shower * 

Philip Stephens and Andre van Hameren 

(N ' 



O 
Oh 



Institute of Nuclear Physics, Polish Academy of Sciences 
ul. Radzikowskiego, 31-342, Krakow, Poland 



(N 1 February 2, 2008 



Abstract 



> 

o 

£f~^ • Presented here is a technique of propagating uncertainties through the parton 

shower by means of an alternate event weight. This technique provides a mechanism 
to systematically quantify the effect of variations of certain components of the parton 
shower leading to a novel approach to probing the physics implemented in a parton 
shower code and understanding its limitations. Further, this approach can be applied 
to a large class of parton shower algorithms and requires no changes to the underlying 



implementation. 

43 ' 



1 Introduction 

As we enter a new era of particle physics, precise knowledge of quantum chromody- 
namics (QCD) will become increasingly important in order to understand the physics 
beyond the standard model. Currently, one of the most useful tools for studying QCD 
is the parton shower approximation. This tool provides a mechanism to connect 
few-parton states to the real world of high-multiplicity hadronic final states while 
retaining the enhanced collinear and soft contributions to all orders. 

Use of parton shower Monte Carlos (MC) has become common-place. Often, 
when one needs an estimate of the uncertainty of a MC prediction several different 
MC programs are used and the differences between them is considered the error pQ. 
Though this technique of estimating the error of the MC is generally acceptable, it 
does little to provide insight into the physics. It has been shown [2] that the un- 
certainties in both the perturbative expansion and the parton distribution functions 
indeed can lead to effects of the order of ten percent. We propose here a technique in 
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which the known uncertainties of the physics can be propagated through the parton 
shower framework. This technique provides alternate weights to an event generated 
by a MC without having to change the basic structure of the MC program. We 
feel this technique could be valuable when determining how various improvements in 
the parton shower will impact the MC predictions. Furthermore, this gives a more 
satisfactory description of the errors in a MC prediction. 

We begin by applying the variational technique to the parton shower probabil- 
ity densities. Using this technique we are able to define the appropriate weights 
associated with the variations. We then define a MC implementation and present 
numerical examples of applying this technique to the parton shower MC. We present 
the algorithm and formula for a variation to the running coupling and the structure 
of the kernel to show the method works as expected. We discuss a way that this pro- 
cedure may be able to be used to estimate the effects of next-to-leading log (NLL) 
terms, additionally we consider how to use this technique to map between two parton 
shower implementations which rely on different interpretations of the kinematics. 



2 Variation of Parton Shower 

In many parton showers [31 IU [6] one starts with the fundamental probability 
density (for one emission) defined as 



V = f R (y)exp[-J d n y'f v (y')\. (1) 



Here the function fn(y) is the distribution of the real emission while fv{y) is the 
virtual contribution. In both cases the precise definition of y is specific to the imple- 
mentation. Furthermore, the limits of integration in the virtual component are also 
specific to the implementation: how the infra-red limit is treated, the definition of 
resolvable versus unresolvable emissions and the ordering of variables. 
For a time-like shower f R = fy and is given by 

us) = ^mptsi. (2) 

where g(y) is some abstract function used to determine the scale of the running 
coupling. We find a similar result for the constrained MC [5] ; for a space-like shower 
using the backward evolution algorithm we find Jr = and 

fv®*)= ^ mn ;/ z '* pw, (3) 

2vr f{x,y) 

where f(x,y) is the PDF at energy fraction x and scale given by some combination 
of the components of y. We can explicitly see that one of the components of y is z, 
a momentum fraction. 

In the forward (time-like) evolution algorithm, as well as the non-Markovian al- 
gorithm, P(y) is just the Alteralli-Parisi [7] splitting function divided by the scale. 
In the numerical results here we consider only the forward evolution algorithm. 
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We now define a functional to represent our functions and fy(y) 

FrMv)} = f R (y) ; F v [<pm = My)- (4) 

Here (f(y) are the functional components of F R jy which we want to vary (e.g. the 
running coupling or the kernel). This defines the distribution of one branching as 

V[<p(y)] = F R [cp(y)) exp I - J d n y'F v [^y')} \ . (5) 

We can find the variation of this by 

SV = V[(cp + 6ip)(g)]-V[<pm- (6) 

If we define 

5F R/V = F R/V [{<p + dcp)(y)} - F R / V [ip(y)\, (7) 

then 



5F R 



r<{y) 

5V = V[l + ^)exp\- I d n y'5F v }-V, (S) 



from which we have a weight 



V + 5V / 5F R \ ( 



w = — ^— = \1 + y^j exp I - J d n y'SF v j . (9) 

If, at first, we assume that F R / V is proportional to a linear product of functions, 
(p a , and we consider variations of only one function we find 

SF aR/v (y) = d a F R/v [<p($)]6ip a (y), (10) 

with 

d a FMy)] = ^M. (11) 

We now turn to varying multiple functions at the same time; again assuming 
that each F R /y is proportional to the linear product of all ip a . At lowest order in 
variations 

thus the weight is given by 



w 



1 + E d f^m 6tPam ) 6XP (" E f d n y'd a F v My')]SMy)) • (13) 



One could keep higher order terms of the variations, but the general formula have 
not been presented here. The weight for such results, given in eqn. Q, would be the 
same with the modified definitions of 5F R /y. 
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Figure 1: Final-state parton branching. The blob represents the hard sub-process. 



The weights defined in eqn. ([9|) are relative to the original probability density for 
one emission. To get the total weight for the full event, we must consider 



and thus 

SV E =V e [<P + Sep, {yi}] - V E [p, {&}]. 
This leads to a total event weight given by 

V E + SV E 



w e 



n 



W: 



(14) 
(15) 
(16) 



3 Example Parton Shower Kinematics 

For the examples given here we will use as a model bremstrahlung emissions from 
one quark line. This is shown in fig. [IJ For the numerical results presented in the 
following sections we now define a concrete implementation of the kinematics of the 
parton shower. We will use the variables of Herwig++ [31 [8]. 

We begin by looking at the ith gluon emission — > qi + ki. In the Sudakov 
basis this is 

qi = aip + (3in + q ±i , (17) 

with p the jet's "parent parton" momentum and n a light-like "backward" 4-vector. 
These obey 

p 2 = m 2 ; n 2 = ; p ■ n = 1 ; q ±i ■ p = q ±i ■ n = 0. (18) 



We find 



Pi 



fi -q±- afm 2 



2ctip ■ n 

Lastly we define the momentum fraction and relative transverse momentum as 

ou 



Oii-l 



-, P±i = Q±i ~ ZiQu-l- 



One then finds 



i 2 



+ 



pi, 



Z{ 1 Z{ ( 1 Zi) 



(19) 
(20) 
(21) 



The evolution variable used for these examples is 



Pi 



Ql 



z 2 (l-z) 2 z 2 z{l-z) 2 ' 

where \x = max(m, Q g ) and Q g is a cutoff parameter of the model. 
Here the splitting kernel, in the quasi-collinear approximation, is 



(22) 



P gq (z,p 2 ± ) = C F 



1 + z 2 2z(l-z)m 2 
1 — z p 2 ± + (1 — z) 2 m 2 

With the definition eqn. (|22|) the branching probability is 

, _ 2l dq 2 dz 



(23) 



<IBU,- <,,;) =^-as[z 2 (l-z) 2 q 2 ] , 



q 2 1 — z 



1 + z' 



2m 2 



zq z 



(24) 



The phase-space boundaries are defined by requiring a real transverse momentum, 
this is found from eqn. (i22|) . We denote the solutions as z^^ 2 ). 

Using the branching probability we can now define the probability density of our 
parton shower. Here we find 



p = dB(q^qg)exp{- J dB{q -»■ qg)j , 



(25) 



which gives 



FrMv)} 



— a s [z (1-z) q ]^ 
F[(a s ,P qq )(z,q 2 )], 



2m 2 



1 + z 2 



1 — z z(l — z)q 



7,2 



1 



2nq' 



as(z,q )P qq (z,q ) 



(26) 



4 Running Coupling Uncertainty 

The first variation we wish to study is that of the running coupling. The running 
coupling has several sources of uncertainty. For example, one may wish to study 
how the variation of the argument to the coupling changes the results of the MC. 
This is often done as a practical way to control the perturbative series of the running 
coupling. 

To illustrate the method described in the previous section, we choose the variation 
of the coupling to be due to the uncertainty in the measurement of the coupling. 
Standard practice in MCs is to take the central value of the coupling at m 2 z and use 
the two-loop renormalization group equation (RGE) to run the coupling to alternate 
scales. At the Landau pole, the perturbative series breaks down and one imposes, 
by hand, some treatment of the coupling below that scale. We propose to trace the 
uncertainty in the value of the running coupling at m 2 z through the RGE and to 
study the effect of this variation on the predictions of the MC. 
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Qualitatively, one would expect to see that if we take the upper bound of the 
running coupling value that we have more emissions and overall a harder p± spec- 
trum of the outgoing quark. If the lower bound is taken, we expect the opposite. 
The magnitude of the fluctuations are governed by the evolution of this uncertainty 
through the RGE. This is exactly what we will see. 

We want to stress that the method described in the previous section is not appli- 
catble only to this particular choice of running coupling variance. One could use any 
choice for the uncertainty of the running coupling and the mechanism of propagating 
this through the parton shower is unchanged, another option is 

5a+(Q 2 ) = a s (Q 2 ) - a s (2Q 2 ) ; Sa^Q 2 ) = a s (Q 2 ) - a s (Q 2 /2), (27) 

for example if one wanted to estimate the error due to the truncation of the pertur- 
bative series for as- 

In most modern applications one uses the two-loop RGE for the running coupling. 
For illustrative purposes in the following examples we choose the running coupling to 
take the value asim 2 ?) = 0.117 ± 10%. We then use the running defined by two-loop 
RGE to run to all scales. In practice we solve for Aq CD which gives the correct value 
of as{m 2 z ) and find the upper and lower variances of Aq CD to give the 10% variation. 
These are then used as input to give the running coupling at alternate scales; mz is 
chosen as 91.2 GeV. 

In order to use the running coupling we must also specify its behaviour below the 
Landau pole (Aqcd)- For these examples the choice made is to freeze the value of 
the coupling at 500 MeV. 

For the two-loop running coupling we have 

4vr 01, f 0ias(Q 2 ) \ a , ( Q 2 



with 



. PI , / W«SW ) \ _ r, , V / oe x 



A) = n-ln f , (29) 



(1"2 v"/)- ( 3 °) 



3 

3 

where n/ is the number of flavours; for these examples this is fixed at 5. 

This can be solved numerically both for the choice of Aq CD and for the value of 
the running coupling. In this case we find Aqcd to be 230^^1 MeV. Of course, there 
are many technical issues related to the running coupling. For example matching the 
running at the threshold of heavy quark masses. For simplicity we ignore these issues 
in the numerical results that proceed. 

We now define the variation in the running coupling to be 

5a+(Q 2 ) = as(Q 2 ;A 2 QCD + SA 2+ )-as(Q 2 ;A 2 Q c D ) (31) 
5a s {Q 2 ) = a s (Q 2 ;A 2 QCD -5A 2 -)-as(Q 2 ;A 2 QCD ). (32) 
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Figure 2: The value of the running coupling and the bounds given by the input at m 2 z . 

With this we can define the weight for each emission due to the variation of the 
running coupling as 

Wa < = 1+ 6XP ~ /, ~W \ dz5a s (z,q)P qq {z,q)\ . (33) 

\ V-K^iiHi) / \ Jqt H Jz i (q 2 ) J 

This weight is normalized to the unvaried weight 1 MC event. 

Figure [2] shows the range of values for the running coupling for the two- loop RGE 
with a 10% variation in the value at m z . We can see that as one approaches A-q CD 
the variation gets larger. This region with a small scale is often where the emissions 
from the parton shower reside; therefore the uncertainty in the measurement of the 
running coupling can have a non-neglible effect on the MC predictions. 

Figure [3] shows the effect of the variations on a) the number of emissions and b) 
the p^_-spectrum. Each of these plots is divided into two panels. The top panel shows 
the results while the bottom panel shows the ratio of the reweighted vs. unweighted 
MC. Here we see the behaviour one would expect. The higher bound of the coupling 
gives more emissions, and a smaller p^-spectrum while the lower bound does the 
opposite. In this figure we chose a massless quark with initial scale q 2 = (1 TeV) 2 . 

In figure [3b, we have also included in the ratio panel the ratio of the reweighted MC 
vs. an alternate MC sample generated by changing the central value of the running 
coupling to the upper or lower bound of the variance. We see from these results that 
the reweighted shower does produce the same results as reimplementing the shower 
with the changed running coupling. 

5 Kernel Variations 

What may be of more interest from a theoretical side is to vary the structure of 
the splitting kernel itself. For example, one could start with the collinear splitting 
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Figure 3: a) The distribution of the number of emissions for m 2 = and q 2 = (100 GeV) 2 . 
b) The p 2 ± distribution. In both cases the solid line is the unvaried case while the dashed 
in the upper bound and the dot-dashed is the lower bound. Additionally, both figures 
contain a second panel which shows the ratio of the varied to unvaried results. In a) we 
have also included the ratio of the reweighted MC vs. an alternate MC sample generated 
by changing the central value of the running coupling to the upper or lower bound of the 



variance. 
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,%2\ 



kernels and vary them by the mass dependent quasi-collinear kernels to see whether 
such changes introduce dramatic effects on a set of observables. 

The benefit to the procedure presented here is that there is no need to change the 
fundamental structure of a given MC. In fact one could add an option to their code to 
keep track of the alternate weights, without changing at all their basic MC program 
logics and structures. One caveat is that though this method will give an accurate 
estimate of the variations given, this is only true for regions of phase space in which 
the original MC fills. If some regions of phase space are empty, or rarely entered, the 
changes in that region due to the variation will still lack significant statistics. 

For completeness, we present a study which shows the effect of introducing the 
quasi-collinear kernel as a variation on the collinear kernel. The collinear kernel is 
simply 

p«W = yzT- ( 34 ) 

To obtain the quasi-collinear kernel, eqn. (|23p . we must define a variance of 

SP qq (z,f) = - m (35) 
z(l — z)q A 

With this variance we find the alternate weight, for the zth emission, is given by 

/ bP qq {z % A?) \ ( f^df ft ^ 
w Pi = ( 1 + ^¥T exp / — / dz a s [z (I - z) q \ 5P qq {z, q 

V Pqq\Z-i,qi ) J \ Jqf Q J z i 

(36) 

and the total weight due to the kernel variation is the product of the weight for each 
emission. This weight is normalized to a weight 1 event with no variations. 

We now show the result of this variation in the fig. 0] when showering a top quark 
with mass 175 GeV from an initial scale of 1 Tev. In fig. we show the effect 
that the quasi-collinear variation has on the distribution of the number of emissions. 
As would be expected, for larger masses we have fewer emissions. Fig. 0b shows 
the p\ spectrum of the outgoing quark. As before these figures are divided into two 
panels. The top panel shows the results while the bottom panel shows the ratio of the 
reweighted MC vs. the unweighted one. Again, in fig. [Jh. the ratio panel also includes 
the ratio of the reweighted MC vs. an alternate MC sample created by changing the 
kernel in the MC to the quasi-collinear kernel. Again we see that this ratio is 1 with 
small variations. 

Though we do not provide any numerical results of varying both the kernel and the 
running coupling simultaneously, we will present the formulae for these weights. In 
the time-like evolution that we have discussed, we can define the weight for emission 
i from eqn. (fT3|) as 

± _ /, 6ajj(zi,tfl) , 5P qq (zi,qfy\ 

)) 

dz (P gq (-)6a±(-) + 6P qq (-)a s (-)) 



as(zi,tf) P qq (zi,qf) J 



x exp 



^dq 2 



Q 2 Jz-tf) 
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Figure 4: a) The distribution of the number of emissions for the collinear kernel and the 
quasi-collinear kernel for m 2 = (175 GeV) 2 and q 2 = (ITeV) 2 . b) The distribution of 
the pi of the outgoing quark for the collinear and quasi-collinear cases under the same 
conditions as a). The solid line shows the result when the quasi-collinear kernel is used, 
the dashed line shows the result when the variation in eqn. fl35l) is applied and the events 
are weighted. Again, the second panel shows the ratio of the varied to the unvaried MC. 
In figure a) also included is the ratio of the reweighted MC vs. an alternate MC sample 
created by changing the kernel in the MC to the quasi-collinear kernel. 



10 



and the weight for the event by we = Yli w f- I n e 1 n - the • represents the pair 



5.1 Combining Kernel with Running Coupling 

Another potential variation that may be of interest is to vary the kernel by a term 
proportional to the running coupling. Such a variation could be used to introduce 
some NLO effects into the kernel. 

We must now determine the appropriate variations in this case. We start with 
the general case for the sum of terms up to a$ 

F[(a s ,pM,...,pW)(z,?)] = £ r^p-) P^W), (38) 

i=i ^ ' 

of which we now have N + 1 functions to vary over. We compute 5F[np] = F[ap + 
Sip] — F[np] and keep the lowest order in the variations. 



8F[tp] = ^2(a s + Sa s y (P® + 6P^J - JV s P (i) , 

i i 

~ Yl ( iai s l6a s p(i) + • (39) 



which is equivalent to the functional derivative. If we keep higher order terms the 
equation is 



N 



sf=y: 



i=i 



(40) 



where N is the length of the vector <p. 

For the examples given here we will look at only N = 2 and set P^ = and 
study the variation around that choice. This means we have 

5F = 5a s + £PW) + [2a s Sa s + 5a 2 s ] 6P$ + a s 6P$ + a 2 s 5P$ . (41) 

If we consider only the lowest order in the variations 

6F « 6a 8 P$ + a s 6P$ + a 2 s 6P$ . (42) 

(2) 

We now turn to an example. We choose the form of SPqq (z) according to full NLO 
kernel [9l[T0]. This is composed of two parts, the flavour singlet (S) and non-singlet 
(V) contributions 

SPg\*,?) = PS?Hz) + P!FH*), (43) 

where these functions are defined in the appendix. We choose 5P^ = and 5a = 
for these examples. 
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Figure [5^, shows the effect on the number of emissions and fig. [5b shows the effect 
on the p^-spectrum of the outgoing quark line. We see that the number of emissions 
is slightly higher with a harder spectrum. 

The construction of a next-to-leading log (NLL) parton shower has the problem 
of negative values for the splitting kernels. These destroy the probabilistic interpre- 
tation of the Sudakov form factors. Naively, one would assume that this will destroy 
any meaningful results for the NLL weights. In our case, this is not true. We are 
reweighting the total density according to the NLL corrections. These may intro- 
duce large or negative weights to the reweighted shower, but this is necessary as this 
correctly describes the density. In the inclusive picture, these negative weights are 
integrated over and pose no problem; exclusively, these negative weights must be 
treated correctly in the analysis. 



6 Variation of Kinematics 

We now consider another use of the alternate weights. Here we wish to use these 
weights to transform one parton shower into another. This, of course, is not an 
exact transformation. This requires additional knowledge about the structure of the 
alternate parton shower. 

The idea is to use the variables generated by one shower and reshape the distri- 
bution to give the results if an alternate shower was used. In this section we discuss 
the intrinsic kinematical definitions. 

Consider a new kinematics, similar to the one used in Pythia [3]. Here we wish 
to order the parton shower in virtuality (Q 2 ). This requires a mapping from q 2 into 
Q 2 . First we have the definition of the transverse momentum as 

pl = z(l-z)(Q 2 -Q 2 ), (44) 



where z is the momentum fraction in the Pythia-like kinematics. From eqn. (|22 j) we 
find (neglecting Qq) 

q2 = * 2 (i-*)V-M 2 (i-*) 2 -*Qg (45) 

z(l - z) 

We also have the boundary of real phase space given by the requirement that there 
is a real ~p\ and the imposition of a particular ordering scheme. 

There is a different interpretation of the meaning of the momentum fraction z in 
the Pythia-like and Herwig-like shower; they have the same distribution, however. 
We compensate for this by constructing the full four-momentum from the Herwig-like 
shower and deconstructing the associated variables for each emission. The weights 
can then be computed from this. This method has the additional benefit that the 
four momentum configuration is identical in both cases; thus hadronization effects 
and hadron decays are identical. 

We now turn to the structure of the probability density itself. For our original 
kinematics we find (for the massless case) 

F[(a s ,P qq )(z,q 2 )} = *s(z 2 (l -z) 2 q 2 )P q<i (z) , _ f _ _ 

Z7T q z 
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Figure 5: a) The distribution of the number of emissions using the collinear kernel at 
O(ots) and applying the variation discusses in the text at 0{a 2 s ). b) The pi distribution 
of the outgoing quark under the same conditions as a). 
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while in the Pythia-like shower we have 
F[(a s , P qq )(z, Q 2 )] = ~ JV ~ v ; '-^^ 0{QU ~ Q Z )0{z + - z)0{z - f"). (47) 



a s (z( l - z)Q 2 ) P qq (z ) Q ^ 2 
2tt Q 2 



From these we can define our variations such that 

F[(a s , P qq )(z, Q 2 )} = F[(a s , P qq )(z, q 2 )\ + 5F, (48) 

thus 

P(z, Q 2 )dQ 2 dz =(V + SV)(z, q 2 )dq 2 dz. (49) 

From this we find 

SF = F[(a s , P qq )(T(z, q 2 ))]J(z, Q 2 ) - F[(a s , P qq )(z, q 2 )}, (50) 

where J is the Jacobian factor for the coordinate transformation T(z, q 2 ). These are 
defined in the appendix. At this point we can exploit the analytic structure of the 
Sudakov form factor, 

A(t;*o) = A(*;ti)A(ti;*o). (51) 

This allows use to seperate the weights into the real and the Sudakov components 
and to calculate the Sudakov components over the full evolution scale, rather than 
just the scales between each emission. This gives 



^ " a„«l,,~§) ■ (52) 

The resulting weight for the real emissions is 

jr) _ a s m - -zdQ 2 )P^\zi)0{QU - QMz + - - n g 2 J& Q 2 ) (53) 

a s {z 2 {l- Zl ) 2 q 2 )P^\z,q 2 ) ~~ ' Q 2 

Here the 9 functions for the Herwig-like evolution are ignored as they are always 
fulfilled by the original shower construction. The total weight is given simply as 

N 

The question now is what does the weighted shower physically give us? This 
gives us the weight, relative to the unweighted original shower, of producing the 
kinematical configuration via the other shower. For our example here this means 
that it will weight our Herwig-like shower to be that of the Pythia-like construction. 
Our weighted shower will produce events that are both ordered in virtuality and in 
angle. Comparing the weighted results versus an independent implementation of the 
full Pythia-like shower would illustrate, for any observable, the effect of the different 
limits in phase-space inherent in each implementation. Furthermore, it could be 
used to illustrate the effects of alternate choices of ordering; e.g. colour connections 
between jets. 
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Figure 6: 1 — T for the Herwig-like shower and reweighted to a Pythia-like shower, as 
described in the text. These differences are due to the different kinematics definitions used 
in each shower. The bottom panel shows the ratio of the Pythia-like vs. Herwig-like. 

To illustrate this technique we use as a model e + e— annihilation into a qq pair. 
This pair then undergoes final state radiation, but the subsequent emissions do not. 
We reconstruct the kinematics of the event and, in order to conserve y/s, we rescale 
each jet by a common factor, k, such that 



N 



Qi+ k Ph 



(55) 



where qf is the virtuality of jet i. To illustrate the reweighting between the Herwig- 
like and Pythia-like shower we study the thrust observable. This is given by 



n 



max • 

n 



Ef=i \Pi 



(56) 



This observable was chosen as the thrust has a strong correlation to the hardest 
emission, but also is effected by subsequent emissions. As we do not shower the 
emitted gluons, studying an observable which have a strong dependence on 2 or 
more emissions is not as illustrative. 

Figure [6] shows the result for = 1 TeV. We can see the deviations, and as 
expected they are not too large. As these are not the result of a full event generation 
it is not useful to compare these to data. 

7 Conclusion 



We have presented a new approach to understanding the errors associated with a MC 
prediction. This approach can be added to almost all currently existing MC programs 
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without changing the physics or the behaviour of the code. Instead, we have provided 
a method to track alternate weights for events. These alternate weights provide the 
tool to reshape MC predictions to see what such a prediction would be if various 
pieces of the MC were altered. 

Though this technique is quite successful, it cannot compensate for all possible 
alterations. As this algorithm provides an alternate weight for an event generated 
by a MC it cannot provide events which cannot be generated by the original MC. 
This means that some of the physical limitations of an already existing code cannot 
be overcome through this method. We do not see this as a drawback, however. The 
purpose of this technique is to understand the physics and the limitations inherent 
in a MC implementation. To this end, such limitations of this technique can provide 
valuable insight. 

This paper has provided numerical examples of a toy parton shower model based 
on the real MC behaviour of Herwig++ It may be quite illustrative to apply this 

method to a fully featured general purpose MC, including hadronization and hadron 
decay, to see how much variation exsists in such a parton shower implementation. 
With such an implementation one may be able to check the accuracy of many MC 
predictions and to understand the limitations of these predictions. 
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A NLO splitting function 

Here we present the formulae for the NLO splitting functions used in this paper. 
These are defined in the MS factorization/renormalization scheme. First we present 
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the flavour singlet contribution [9] 

P^ 2) {z) = C F l-l + z + ^ [(1-3^) In 2- {l + z)\n 2 z\ 



— In z + 2 In z ln(l — z) 



Pqq(z) + 2p qq (-z)S 2 {z) 



14 

+ ('/■■(' y< y (!-«) + 



11 



1 2 

In z H — In z + 



67 
18 



7T 

~6~ 



P<ra(*) 



-p OT (-*)fl2(*)} 

CFT f {-% + f.+ (l0z + f, 2 + 2)ln, 



112 o 40 , „ 
-— . 2 + --2(l + z)ln 2 z 



10 2 , 

1 — In z 

9 3 



In this formula we have 

Pqq{ Z ) = 

S 2 (z) = 



C F 



1 - z 



1-z, 



i+z dz , ( 1 — 
— m 

z \ z 

l + z X 

-2Li 2 (-z) + - In 2 2 - 2 In z ln(l + z) 



IT 



Ca = 3; T 



f 



Now we give the flavour non-singlet function [TO 



pV(2) 

<?<? 



Cl 



2 In z ln(l — z) + — In z 



Pqq(z) 



+ lnz-i(l + z)ln 2 2-5(l 



CfCa 



1 



11, 

In z H In z + 



6 



67_^ 
18 6 



P<ra(*) 



20 



(1 +z)lnz + — (1 - z) 



+C F T f 



2 , 10 

— In z H 

3 9 



Pqq( Z ) ~ ^ C 1 - *)} • 



(57) 

(58) 



(59) 
(60) 



(61) 



We note that the superscripts given here differ from the normal convention. Here 
they indicate the total number of powers of as in the branching probability. Normal 
convention decrements these by one to indicate the total order of expansion in the 
splitting kernel. 

B Coordinate Transformation Jacobian 

In this appendix we compute the full Jacobian factor for the transformation between 
the evolution variables used in the Herwig-like shower and those used in the Pythia- 



17 



like shower. This transformation has the form 

T(q 2 ,z)^(Q 2 ,z). 



(62) 



This is to be done after the full momentum reconstruction so we know all of the 
components of the momenta. We need to numerically evaluate the Jacobian factor 
for the weights of the real emissions. 
We start with 

ffi a 2 O 2 



and 



Together we find 



We now find z as 



Q 2 i = ill 



z 2 (l-z) 2 z 2 z(l-z) 2 ' 

k 2 



z 



A 



l-z z(l-z)' 

2 _ Qi+(z-l)fi 2 2 _ 



Q = 



+ q 2 z(l - z). 



q° ctip + fan zai-ip° + fan 



„0 



„0 

Qi-1 



„0 

Qi-i 



From the Herwig++ variables we have 

.,2 , 2 

fa 



..2^,2 



2atip • n 



Qt ~ Ql 



2 2 

z ai_ x m 



We now use 



2zctj_ip • n 



q±i = p±i-zq±i-i, 

7 = 2ai^p ■ n, 



to find 



fa = —[q 2 +q 2 z 2 (l-z) 2 -v 2 (l-z) 2 + zQ 



+2z^q 2 z 2 (l - z) 2 - fi 2 (l - zf + zQ 2 \q ±i -i\cos 9 



, -2 2 2 21 

+q±i-! - ZiOti-xm \ 



Using these formulae we are now able to compute the full Jacobian 

J(z,Q 2 



^ = d(Q 2 ,z) 
d(q 2 ,z) 



dQ 2 dz 

dq 2 dq 2 

dQ 2 dz 

dz dz 



(63) 
(64) 
(65) 

(66) 

(67) 
(68) 



(69) 
(70) 



(71) 



(72) 
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We have 

dQ 2 



dq 



2 



z(l - z), (73) 



g - 2^ + ,74, 

dz Pi n° 



dq 2 dq 2 g?_ 1 ' 



(75) 



dz ai_ip° d@i n° 

Tz = -^ + — — • 

9 Pi 2/, n2 / ■, , COS 




9 « 2 V ^ 2 ^ 2 (1 - ^) 2 - M 2 (l - z) 2 + zQ2 

9 2 (1 - A* 2 (z 2 - 1 



<9z 7Z 2 7 7 \ z 2 



(77) 
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g 2 2z(l-z)(l-2z) + 2^ 2 (l-z)+Q 9 a^m 2 



"f^q 2 z 2 (l - z) 2 - n 2 (l - z) 2 + zQ 



|g_u_i|cos( 
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We can now completely determine the Jacobian factor for transforming between 
the two kinematics. Thus the weight for the real emission is 



w (R) = as(z(l - z)Q 2 )P£\*, Q 2 )J(z, Q 2 ) Q 2 m 
a s (z 2 (l-z 2 )q 2 )p£\z,q 2 ) ' 



Finally, it is left to compute the cutoff virtuality from the cutoff in q 2 . Since we 
know the outgoing quark has virtuality /j, 2 we assume that this is the cutoff on the 
virtuality ordered shower. One could find from 

the value of z to use in 

Q 2 = q 2 z(l-z)+fi 2 . (81) 
The weight for the Sudakov factor is then 

P \ hi ^J Z -(Q 2 ) aZ 2^ ^qq ) 

w A = ^ J f . (82) 

6XP l k 2 o ^JzZiq 2 ) ^ m { ,q) 
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